

* PROBABILITY OF GETTING MARRIED 
qui {
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if AGE1<65
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
drop if labor_force1==.

drop if year==2013 // we don't want to use data from 2014 which is the second year of the panel started in 2013


gen Mar_Indic1 = 0 if MAR1==0 & MAR3==0
replace  Mar_Indic1 = 1 if MAR1==0 & MAR3==1

	
	label define incquint 1 "Income Quintile=1st" 2  "Income Quintile=2nd" 3  "Income Quintile=3rd" 4  "Income Quintile=4th" 5  "Income Quintile=5th"
label values inc_quintile_Y1  incquint
label var Mar_Indic1 "Single to Married"

eststo clear
	
replace labor_force1=1 if (labor_force1==2 | labor_force1==3)

logit Mar_Indic1 AGE1 AGE1_sq AGE1_cub     i.inc_quintile_Y1   i.labor_force1   // 
eststo
	esttab using "${out_tables}/Reg1.tex", varwidth(25) nogaps   compress label nobaselevels indicate( "Cubic Age = *AGE1* *AGE1_sq* *AGE1_cub*" ) replace  ///
	 nonumbers noconstant eqlabel(none) star(* 0.1 ** 0.05 *** 0.01) ///
	se   b(3) 
	eststo clear

	fillin AGE1 inc_quintile_Y1 labor_force1
	replace AGE1_sq=AGE1^2 if _fillin==1
	replace AGE1_cub=AGE1^3 if _fillin==1
	
	* save predictions
logit Mar_Indic1 AGE1 AGE1_sq AGE1_cub     i.inc_quintile_Y1   i.labor_force1  
	predict Y, pr
	collapse Y, by(AGE1 labor_force1 inc_quintile_Y1)
	* add all the "o" states
	 reshape wide Y, i(AGE1 inc_quintile_Y1 ) j( labor_force1)
	 gen Y2=Y1
	 gen Y3=Y1
	 reshape long Y, i(AGE1 inc_quintile_Y1 ) j( labor_force1)
	sort AGE1 labor_force1 inc_quintile_Y1
	drop if AGE1==24
drop AGE1 labor_force1 inc_quintile_Y1

outsheet using "${data_model}\Marriage1.txt", nolabel nonames replace

}


* INITIAL DISTRIBUTION for BEING MARRIED AT 25, By Education, FOR MODEL
qui {
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.

keep if AGE1==25 | AGE1==26
drop if H1==.
table educ_group H1 [aw=Sample_Weight] , stat(mean MAR1)
collapse (mean) MAR1 [aw=Sample_Weight], by(educ_group H1)
sort educ_group H1
drop educ_group H1

outsheet using "${data_model}\Marriage0.txt", nolabel nonames replace
}


* probabilities of separating for model
qui {
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if AGE1<65 & AGE1>=25
drop if labor_force1==.
drop if H1==.
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
drop if year==2013 // we don't want to use data from 2014 which is the second year of the panel started in 2013

gen Mar_Indic2 = 0 if MAR1==1 & MAR3==1
replace  Mar_Indic2 = 1 if MAR1==1 & MAR3==0

	label define incquint 1 "Income Quintile=1st" 2  "Income Quintile=2nd" 3  "Income Quintile=3rd" 4  "Income Quintile=4th" 5  "Income Quintile=5th"
label values inc_quintile_Y1  incquint
label var AGE1 "Age"
	
	label var Mar_Indic2 "Married to Single"
	eststo clear

replace labor_force1=1 if (labor_force1==2 | labor_force1==3)

logit Mar_Indic2 AGE1    i.H1 i.inc_quintile_Y1  i.educ_group  i.labor_force1  // 
eststo
	esttab using "${out_tables}/Reg2.tex", varwidth(25) nogaps   compress label nobaselevels  replace  ///
	 nonumbers eqlabel(none) noconstant  star(* 0.1 ** 0.05 *** 0.01) ///
	se   b(3) 
	eststo clear
		
	fillin AGE1 H1 inc_quintile_Y1 educ_group labor_force1
	
logit Mar_Indic2 AGE1    i.H1 i.inc_quintile_Y1  i.educ_group  i.labor_force1  // [pw=Sample_Weight]
	predict Y, pr
	collapse Y, by(AGE1 H1 labor_force1 inc_quintile_Y1 educ_group)
	 reshape wide Y, i(AGE1 H1  inc_quintile_Y1 educ_group) j( labor_force1)
	 gen Y2=Y1
	 gen Y3=Y1
	 reshape long Y, i(AGE1 H1  inc_quintile_Y1 educ_group) j( labor_force1)

	sort AGE1 H1 labor_force1 inc_quintile_Y1 educ_group
drop if AGE1==24
drop AGE1 H1 labor_force1 inc_quintile_Y1 educ_group

outsheet using "${data_model}\Marriage2.txt", nolabel nonames replace
}


* PROBABILITY THE SPOUSE WORKS, CONDITIONAL ON MARRIED. 
* this depends on age, educ, and H. from ages 25-64.
qui {
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
	
keep if AGE1<65 & AGE1>=25
drop if H1==.
keep if MAR1==1 & MAR2==1 

*EMPST
* 1 Employed at RD [.] Interview Date
* 2 Job to Return to at RD [.] Interview Date
* 3 Job During RD [.] Ref Period
* 4 Not Employed During RD [.]

* spouse status
gen Spouse_status = 1 if MAR1==0 & MAR2==0
replace Spouse_status = 2 if MAR1==1 & MAR2==1 & spouse_LF_Y1==0
replace Spouse_status = 3 if MAR1==1 & MAR2==1 & spouse_LF_Y1==1
label define sp_stat 1 "Single" 2 "Married, Spouse N.E." 3 "Married, Spouse E."
label values Spouse_status sp_stat
label var Spouse_status "Marital and Spouse Status"

	
	fillin AGE1 H1 educ_group inc_tercile_Y1 
	replace AGE1_sq=AGE1^2 if _fillin==1
	replace AGE1_cub=AGE1^3 if _fillin==1


	label var AGE1 "Age"
	label var AGE1_sq "Age sq"
	label var AGE1_cub "Age cub"
	
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==1 
eststo
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==2
eststo
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==3 
eststo
	esttab using "${out_tables}/Spouse_emp.tex", varwidth(25) nogaps   compress label nobaselevels  replace  ///
	nodepvars nonumbers  mtitles("(HS or Less)" "(S College)" "(College)" ) star(* 0.1 ** 0.05 *** 0.01) ///
	se   b(3) 
	eststo clear
	
	
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==1 
predict Y1 if educ_group==1 , pr
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==2
predict Y2 if educ_group==2 , pr
logit spouse_LF_Y1 AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  if educ_group==3 
predict Y3 if educ_group==3 , pr
egen Y=rowmax(Y1 Y2 Y3)
	collapse Y, by(AGE1 H1 educ_group inc_tercile_Y1 )
	sort AGE1 educ_group H1 inc_tercile_Y1 


drop AGE1 educ_group H1 inc_tercile_Y1 

outsheet using "${data_model}\Prob_spouse_work.txt", nolabel nonames replace	

}


* SPOUSE'S INCOME
* for working ages, depends on age, educ, and if she works. 
qui{
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
keep if AGE1<65 & AGE1>=25
drop if H1==.
keep if MAR1==1 & MAR2==1 

drop if spouse_LF_Y1==.

	// drop the extreme values for spouse_WAGE_YR1
	// drop top 1%
	tabstat spouse_WAGE_YR1 if spouse_LF_Y1==0, stat(min p25 p50 p90 p95  p99 max)
	drop if spouse_WAGE_YR1>113.89 & spouse_LF_Y1==0 // top 1%

	fillin AGE1 H1 educ_group spouse_LF_Y1 inc_tercile_Y1
	replace AGE1_sq=AGE1^2 if _fillin==1
	replace AGE1_cub=AGE1^3 if _fillin==1
	
	label var AGE1 "Age"
	label var AGE1_sq "Age sq"
	label var AGE1_cub "Age cub"
	
	
	reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight]  if spouse_LF_Y1==0 & educ_group==1 //[pw=Sample_Weight]
	eststo
	reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1  i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==0 & educ_group==2 //[pw=Sample_Weight]
	eststo
	reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  [pw=Sample_Weight] if spouse_LF_Y1==0 & educ_group==3 //[pw=Sample_Weight]
	eststo
	qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==1 //[pw=Sample_Weight]
	eststo
		qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1  [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==2 //[pw=Sample_Weight]
	eststo
		qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==3 //[pw=Sample_Weight]
	eststo
	esttab using "${out_tables}/Spouse_income.tex", varwidth(25) nogaps   compress label nobaselevels  replace  ///
	nodepvars nonumbers noconstant  mtitles("NE, <=HS" "NE, SC" "NE, Col"  "Emp, <=HS" "Emp, SC" "Emp, Col" ) star(* 0.1 ** 0.05 *** 0.01) ///
	se   b(3) 
	eststo clear

reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==0 & educ_group==1 
predict Y1 if spouse_LF_Y1==0 & educ_group==1 
reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1  i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==0 & educ_group==2 
predict Y2 if spouse_LF_Y1==0 & educ_group==2 
reg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1  i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==0 & educ_group==3 
predict Y3 if spouse_LF_Y1==0 & educ_group==3 


qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==1 
predict Y4 if spouse_LF_Y1==1 & educ_group==1 
qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==2
predict Y5 if spouse_LF_Y1==1 & educ_group==2 
qreg  spouse_WAGE_YR1  AGE1 AGE1_sq AGE1_cub i.H1 i.inc_tercile_Y1 [pw=Sample_Weight] if spouse_LF_Y1==1 & educ_group==3 
predict Y6 if spouse_LF_Y1==1 & educ_group==3 


egen Y=rowmax(Y1 Y2 Y3 Y4 Y5 Y6)

	collapse Y, by(AGE1 H1 educ_group spouse_LF_Y1 inc_tercile_Y1)
	replace Y=0 if Y<0
	replace Y=(Y*1000) /5200 // convert to model units
	sort AGE1 educ_group H1 spouse_LF_Y1 inc_tercile_Y1
drop AGE1 educ_group H1 spouse_LF_Y1 inc_tercile_Y1

outsheet using "${data_model}\Income_spouse1.txt", nolabel nonames replace	


}


* SPOUSE'S MEDICAL COSTS working ages
qui{
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if AGE1<65 & AGE1>=25
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
keep if MAR1==1 & MAR2==1 // & MAR3==1

* HI = health insurance status 

gen spouse_HI = 1 if spouse_OFFER1==1 & spouse_PRIVAT1==1 // had her own offer and is insured
replace spouse_HI = 2 if spouse_PRIVAT1==0 & spouse_MCAID1==0 // not insured
replace spouse_HI = 3 if spouse_OFFER1==0 & spouse_PRIVAT1==1 // did not have offer, but has PHI from somewhere else
label define ins5 1 "Own ESHI offer" 2 "Uninsured" 3 "Has PHI but not own ESHI offer"
label values spouse_HI ins5

reg spouse_TOTEXP_YR1 AGE1 AGE1_sq i.educ_group i.spouse_HI
replace spouse_HI=1 if spouse_HI==3
reg spouse_TOTEXP_YR1 AGE1 AGE1_sq i.educ_group i.spouse_HI

* justify using round 1 - we get higher charges for "uninsured" 
* so the reasoning is that the spouse may have qualified for Medicaid together with the husband because of her medical expenses later in the year. 
* so we use her insurance status at the beginning of the year.  

drop if spouse_HI==.

gen spouse_charges_use = (spouse_TOTTCH_YR1 + spouse_RXEXP_YR1)*.6 
table spouse_HI, stat(mean spouse_TOTEXP_YR1 spouse_TOTTCH_YR1 spouse_charges_use spouse_TOTSLF_YR1)

drop if spouse_HI==.
	fillin AGE1 educ_group spouse_HI	
	replace AGE1_sq=AGE1^2 if _fillin==1
	
	
//reg spouse_TOTEXP_YR1 AGE1 AGE1_sq i.educ_group i.spouse_HI
  reg spouse_TOTSLF_YR1 AGE1 AGE1_sq i.educ_group if spouse_HI==1 // has PHI
predict OOP_wife if spouse_HI==1

  reg spouse_charges_use AGE1 AGE1_sq i.educ_group if spouse_HI==2 // no PHI
predict OOP_wife1 if spouse_HI==2
replace OOP_wife=OOP_wife1 if spouse_HI==2
drop OOP_wife1
	
	collapse OOP_wife, by(AGE1 spouse_HI educ_group)
	sort educ_group  AGE1  spouse_HI
	preserve 
	* for table with some sample OOP that goes into appendix
	keep if AGE1==25 | AGE1==35 | AGE1==45 | AGE1==55 | AGE1==64
	reshape wide OOP_wife, i(AGE1 educ_group) j(spouse_HI)
	sort educ_group AGE1 // TABLE IN APPENDIX
	restore
drop educ_group AGE1 spouse_HI
replace OOP_wife=OOP_wife/5.2 // convert to model units

outsheet using "${data_model}\OOP_spouse.txt", nolabel nonames replace	
}

	
** spouse's OOP when 65+ 
* everyone has Medicare, so here they just depend on age and educ	
qui{	
clear
use "${data}/MEPS_merged_Health_Whites.dta"
set more off
keep if AGE1>=65 
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
keep if MAR1==1 & MAR2==1 // & MAR3==1

tabstat spouse_TOTEXP_YR1 if AGE1>=65 & AGE1<70
tabstat spouse_TOTEXP_YR1 if AGE1>=70 & AGE1<75
tabstat spouse_TOTEXP_YR1 if AGE1>=75 & AGE1<80
tabstat spouse_TOTEXP_YR1 if AGE1>=80 & AGE1<84


tabstat spouse_TOTSLF_YR1 if AGE1>=65 & AGE1<70
tabstat spouse_TOTSLF_YR1 if AGE1>=70 & AGE1<75
tabstat spouse_TOTSLF_YR1 if AGE1>=75 & AGE1<80
tabstat spouse_TOTSLF_YR1 if AGE1>=80 & AGE1<84


reg spouse_TOTEXP_YR1 AGE1 
eststo
	esttab using "${out_tables}/Spouse_MTC.tex", varwidth(25) nogaps   compress label nobaselevels  replace  ///
	nodepvars nonumbers noconstant mtitles("(1)" "(2)"  "(3)") star(* 0.1 ** 0.05 *** 0.01) ///
	se   b(3) 
	eststo clear

*fillin AGE1	
reg spouse_TOTEXP_YR1 AGE1 
predict OOP_wife
 
	collapse OOP_wife, by(AGE1)
	sort  AGE1  
drop  AGE1 
replace OOP_wife=OOP_wife/5.2 // convert to model units

outsheet using "${data_model}\OOP_spouse_old.txt", nolabel nonames replace	
}


** tables 89 - 91 in the appendix with fraction of spouses working, mean spousal income if working and not working
qui{
clear
use "${data}/MEPS_merged_Health_Whites.dta"

keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.
keep if AGE1<65 & AGE1>=25
drop if H1==.
keep if MAR1==1 & MAR2==1 



* spouse status
gen Spouse_status = 1 if MAR1==0 & MAR2==0
replace Spouse_status = 2 if MAR1==1 & MAR2==1 & spouse_LF_Y1==0
replace Spouse_status = 3 if MAR1==1 & MAR2==1 & spouse_LF_Y1==1
label define sp_stat 1 "Single" 2 "Married, Spouse N.E." 3 "Married, Spouse E."
label values Spouse_status sp_stat
label var Spouse_status "Marital and Spouse Status"

	
table  (educ_group) (inc_tercile_Y1) ,   statistic(mean spouse_LF_Y1)   nformat(%5.2f)  
collect title "Fraction of Spouses Working, MEPS"
collect export "${out_tables}/Spouses_work.tex", tableonly replace

	// drop the extreme values for spouse_WAGE_YR1
	// drop top 1%
	tabstat spouse_WAGE_YR1 if spouse_LF_Y1==0, stat(min p25 p50 p90 p95  p99 max)
	drop if spouse_WAGE_YR1>113 & spouse_LF_Y1==0
		tabstat spouse_WAGE_YR1 if spouse_LF_Y1==1, stat(min p25 p50 p90 p95  p99 max)
	drop if spouse_WAGE_YR1>179 & spouse_LF_Y1==1
	
table  (educ_group) (inc_tercile_Y1) if spouse_LF_Y1==0 ,   statistic(mean spouse_WAGE_YR1)   nformat(%5.2f)  
collect title "Mean Spousal Income if Not Working, MEPS"
collect export "${out_tables}/Spouses_inc1.tex", tableonly replace

table  (educ_group) (inc_tercile_Y1) if spouse_LF_Y1==1,   statistic(mean spouse_WAGE_YR1)   nformat(%5.2f)  
collect title "Mean Spousal Income if Working, MEPS"
collect export "${out_tables}/Spouses_inc2.tex", tableonly replace
}


* table 51 in appendix - labor force of husband conditional on spouse's status
qui {
clear
use "${data}/MEPS_merged_Health_Whites.dta"
keep if INSCOP_Y1==1 // this covers R1-R3 which we are using here.

keep if AGE1<=55 & AGE1>=30
drop if H1==.

* spouse status
gen Spouse_status = 1 if MAR1==0 & MAR2==0
replace Spouse_status = 2 if MAR1==1 & MAR2==1 & spouse_LF_Y1==0
replace Spouse_status = 3 if MAR1==1 & MAR2==1 & spouse_LF_Y1==1
label define sp_stat 1 "Single" 2 "Married, Spouse N.E." 3 "Married, Spouse E."
label values Spouse_status sp_stat
label var Spouse_status "Marital and Spouse Status"

* this is for FT
table educ_group Spouse_status ,  statistic(mean EMPLOYMENT_YR1_v1) nototals nformat(%5.2f) 
collect title "Fraction Employed FT, ages 30-55, by Marital Status, MEPS"
collect export "${out_tables}/Mar_LF.tex", tableonly replace

}